Benchmark Refiner

source: https://www.pnas.org/doi/10.1073/pnas.1921046117#sec-1-2

Étant donné un ensemble d'instances d'une famille de ET putativement apparentées, l'outil RepeatModeler Refiner tente de construire un alignement d'amorçage (seed alignement) de haute qualité et d'en déduire une séquence consensus pour la famille.

seed alignement: c'est la technique dominante pour les comparaisons de séquences génomiques à grande échelle. BLASTN en est l'implémentation la plus populaire. Dans BLASTN, les correspondances exactes de w résidus contigus entre les séquences sont identifiées et ensuite étendues en alignements par programmation dynamique (source).

Ce processus est amorcé en commençant par en effectuant d'abord une comparaison complète par paires de toutes les séquences et en sélectionnant, comme consensus initial, l'instance qui s'aligne le mieux sur toutes les autres.
Il est alors possible que la séquence initiale choisie ne s'aligne pas sur une petite partie des séquences d'entrée. Dans ce cas, les séquences non alignées sont maintenues dans un pool de candidats pour être éventuellement incluses dans d'autres itérations de raffinement du consensus.
Les alignements par paires avec la séquence choisie sont combinés dans l'alignement initial et un nouvel outil d'appel de consensus ajusté aux CpG est utilisé pour générer un consensus mis à jour pour la famille.
L'appel au consensus utilisé diffère d'un appel au consensus standard à règle majoritaire de deux manières: il évalue chaque base ancestrale possible ou code IUB de la colonne d'alignement d'origine en utilisant une matrice de substitution neutre, et il recherche une surreprésentation des produits de mutation CpG communs pour identifier correctement l'état ancestral des dimères CpG. La première étape utilise une matrice qui reflète les modèles de substitution d'ADN neutre observés.

import utils
import plotly.express as px

work_dir = "/home/lerat/TE_Aalb/test/cdhit/benchmark_2/"
data = utils.refiner_table(work_dir)

data

	nb_seq-100	cons_length-100	div-100	nb_seq-150	cons_length-150	div-150	nb_seq-200	cons_length-200	div-200	nb_seq-250	...	div-350	nb_seq-400	cons_length-400	div-400	nb_seq-450	cons_length-450	div-450	nb_seq-500	cons_length-500	div-500
cluster
cluster_9	100	168	0.030120	150	164	0.012195	200	164	0.006098	250	...	0.000000	400	164	0.000000	450	164	0.000000	500	164	0.0
cluster_8	100	859	0.000000	150	859	0.001164	200	860	0.001745	250	...	0.000000	400	859	0.000000	450	859	0.000000	500	859	0.0
cluster_7	100	262	0.358634	150	262	0.358634	200	262	0.358634	250	...	0.361059	400	264	0.361059	450	265	0.000000	500	265	0.0
cluster_6	100	208	0.220619	150	215	0.227642	200	149	0.305164	250	...	0.062718	400	297	0.059233	450	277	0.003610	500	277	0.0
cluster_5	100	618	0.003247	150	618	0.003247	200	618	0.003247	250	...	0.003247	400	618	0.003247	450	614	0.000000	500	614	0.0
cluster_4	100	307	0.397260	150	307	0.399751	200	298	0.413098	250	...	0.004032	400	496	0.004032	450	496	0.004032	500	496	0.0
cluster_39	100	922	0.623491	150	927	0.621832	200	927	0.622240	250	...	0.570071	400	3978	0.014478	450	3967	0.000504	500	3965	0.0
cluster_38	100	160	0.369048	150	176	0.011364	200	176	0.011364	250	...	0.005682	400	176	0.005682	450	176	0.005682	500	176	0.0
cluster_37	100	850	0.005869	150	850	0.004695	200	850	0.004695	250	...	0.004695	400	850	0.004695	450	854	0.000000	500	854	0.0
cluster_36	100	317	0.081761	150	318	0.080063	200	318	0.076923	250	...	0.255385	400	318	0.004710	450	319	0.003135	500	319	0.0
cluster_35	100	577	0.012132	150	577	0.005199	200	577	0.001733	250	...	0.005199	400	577	0.001733	450	577	0.000000	500	577	0.0
cluster_34	100	489	0.382413	150	622	0.386139	200	622	0.386139	250	...	0.330203	400	489	0.378323	450	489	0.002045	500	489	0.0
cluster_33	100	1649	0.293329	150	1647	0.289851	200	1642	0.291209	250	...	0.012797	400	1648	0.077392	450	1634	0.000612	500	1634	0.0
cluster_32	100	583	0.017903	150	583	0.017903	200	573	0.349957	250	...	0.000000	400	590	0.000000	450	590	0.000000	500	590	0.0
cluster_31	100	471	0.372709	150	511	0.005871	200	511	0.005871	250	...	0.001957	400	511	0.001957	450	511	0.000000	500	511	0.0
cluster_30	100	799	0.006258	150	799	0.001252	200	799	0.000000	250	...	0.000000	400	799	0.001252	450	799	0.001252	500	799	0.0
cluster_3	100	650	0.326672	150	650	0.326672	200	650	0.326672	250	...	0.000000	400	651	0.000000	450	651	0.000000	500	651	0.0
cluster_29	100	784	0.000000	150	784	0.000000	200	784	0.000000	250	...	0.000000	400	784	0.000000	450	784	0.000000	500	784	0.0
cluster_28	100	688	0.008721	150	688	0.007267	200	688	0.005814	250	...	0.000000	400	688	0.000000	450	688	0.000000	500	688	0.0
cluster_27	100	334	0.024890	150	352	0.309558	200	352	0.309558	250	...	0.002865	400	349	0.000000	450	349	0.002865	500	349	0.0
cluster_26	100	222	0.264108	150	221	0.000000	200	221	0.000000	250	...	0.000000	400	221	0.000000	450	221	0.000000	500	221	0.0
cluster_25	100	308	0.006494	150	308	0.000000	200	308	0.000000	250	...	0.042386	400	330	0.043887	450	308	0.000000	500	308	0.0
cluster_24	100	2382	0.006061	150	2382	0.006061	200	2382	0.004807	250	...	0.000416	400	2403	0.000000	450	2403	0.000000	500	2403	0.0
cluster_23	100	547	0.027223	150	555	0.046847	200	542	0.376481	250	...	0.009042	400	545	0.372727	450	550	0.011765	500	555	0.0
cluster_22	100	152	0.271565	150	139	0.320000	200	139	0.313333	250	...	0.128205	400	162	0.034056	450	161	0.006211	500	161	0.0
cluster_21	100	482	0.010288	150	483	0.009250	200	482	0.010288	250	...	0.376915	400	489	0.376915	450	490	0.000000	500	490	0.0
cluster_20	100	437	0.011390	150	437	0.011390	200	437	0.011390	250	...	0.002268	400	441	0.002268	450	441	0.002268	500	441	0.0
cluster_2	100	883	0.357913	150	881	0.002839	200	880	0.001136	250	...	0.000000	400	880	0.000000	450	880	0.002273	500	880	0.0
cluster_19	100	356	0.407491	150	503	0.330739	200	503	0.330739	250	...	0.409611	400	349	0.407323	450	529	0.005693	500	525	0.0
cluster_18	100	730	0.044400	150	731	0.035762	200	752	0.033312	250	...	0.003846	400	780	0.003207	450	781	0.002564	500	779	0.0
cluster_17	100	342	0.002924	150	342	0.000000	200	342	0.000000	250	...	0.000000	400	342	0.000000	450	342	0.000000	500	342	0.0
cluster_16	100	3858	0.509391	150	3860	0.508812	200	3859	0.508335	250	...	0.522969	400	4075	0.521251	450	4081	0.522514	500	1360	0.0
cluster_15	100	226	0.347826	150	230	0.340517	200	234	0.008547	250	...	0.008547	400	234	0.000000	450	234	0.000000	500	234	0.0
cluster_14	100	126	0.262500	150	134	0.310976	200	134	0.310976	250	...	0.310976	400	129	0.287926	450	129	0.294118	500	194	0.0
cluster_13	100	262	0.019231	150	258	0.000000	200	258	0.000000	250	...	0.000000	400	258	0.000000	450	258	0.000000	500	258	0.0
cluster_12	100	2111	0.010646	150	2111	0.008753	200	2110	0.006152	250	...	0.003779	400	2118	0.003307	450	2116	0.001418	500	2116	0.0
cluster_11	100	633	0.451852	150	623	0.457330	200	586	0.510676	250	...	0.011121	400	1665	0.007514	450	1661	0.010533	500	1662	0.0
cluster_10	100	310	0.247173	150	309	0.000000	200	309	0.000000	250	...	0.000000	400	309	0.000000	450	309	0.000000	500	309	0.0
cluster_1	100	429	0.026713	150	432	0.030093	200	434	0.027714	250	...	0.013889	400	432	0.009259	450	432	0.009259	500	432	0.0
cluster_0	100	431	0.328125	150	431	0.325893	200	431	0.328125	250	...	0.325893	400	431	0.330357	450	465	0.000000	500	465	0.0

40 rows × 27 columns

Sampled sequences

sampled_seq = data[['nb_seq-100', 'nb_seq-150', 'nb_seq-200', 'nb_seq-250', 'nb_seq-300', 'nb_seq-350', \
                       'nb_seq-400', 'nb_seq-450', 'nb_seq-500']]
sampled_seq

	nb_seq-100	nb_seq-150	nb_seq-200	nb_seq-250	nb_seq-300	nb_seq-350	nb_seq-400	nb_seq-450	nb_seq-500
cluster
cluster_9	100	150	200	250	300	350	400	450	500
cluster_8	100	150	200	250	300	350	400	450	500
cluster_7	100	150	200	250	300	350	400	450	500
cluster_6	100	150	200	250	300	350	400	450	500
cluster_5	100	150	200	250	300	350	400	450	500
cluster_4	100	150	200	250	300	350	400	450	500
cluster_39	100	150	200	250	300	350	400	450	500
cluster_38	100	150	200	250	300	350	400	450	500
cluster_37	100	150	200	250	300	350	400	450	500
cluster_36	100	150	200	250	300	350	400	450	500
cluster_35	100	150	200	250	300	350	400	450	500
cluster_34	100	150	200	250	300	350	400	450	500
cluster_33	100	150	200	250	300	350	400	450	500
cluster_32	100	150	200	250	300	350	400	450	500
cluster_31	100	150	200	250	300	350	400	450	500
cluster_30	100	150	200	250	300	350	400	450	500
cluster_3	100	150	200	250	300	350	400	450	500
cluster_29	100	150	200	250	300	350	400	450	500
cluster_28	100	150	200	250	300	350	400	450	500
cluster_27	100	150	200	250	300	350	400	450	500
cluster_26	100	150	200	250	300	350	400	450	500
cluster_25	100	150	200	250	300	350	400	450	500
cluster_24	100	150	200	250	300	350	400	450	500
cluster_23	100	150	200	250	300	350	400	450	500
cluster_22	100	150	200	250	300	350	400	450	500
cluster_21	100	150	200	250	300	350	400	450	500
cluster_20	100	150	200	250	300	350	400	450	500
cluster_2	100	150	200	250	300	350	400	450	500
cluster_19	100	150	200	250	300	350	400	450	500
cluster_18	100	150	200	250	300	350	400	450	500
cluster_17	100	150	200	250	300	350	400	450	500
cluster_16	100	150	200	250	300	350	400	450	500
cluster_15	100	150	200	250	300	350	400	450	500
cluster_14	100	150	200	250	300	350	400	450	500
cluster_13	100	150	200	250	300	350	400	450	500
cluster_12	100	150	200	250	300	350	400	450	500
cluster_11	100	150	200	250	300	350	400	450	500
cluster_10	100	150	200	250	300	350	400	450	500
cluster_1	100	150	200	250	300	350	400	450	500
cluster_0	100	150	200	250	300	350	400	450	500

The consensus length

cons_length = data[['cons_length-100', 'cons_length-150', 'cons_length-200', 'cons_length-250', \
                       'cons_length-300', 'cons_length-350', 'cons_length-400', 'cons_length-450', \
                       'cons_length-500']]
cons_length = cons_length.T
cons_length

cluster	cluster_9	cluster_8	cluster_7	cluster_6	cluster_5	cluster_4	cluster_39	cluster_38	cluster_37	cluster_36	...	cluster_17	cluster_16	cluster_15	cluster_14	cluster_13	cluster_12	cluster_11	cluster_10	cluster_1	cluster_0
cons_length-100	168	859	262	208	618	307	922	160	850	317	...	342	3858	226	126	262	2111	633	310	429	431
cons_length-150	164	859	262	215	618	307	927	176	850	318	...	342	3860	230	134	258	2111	623	309	432	431
cons_length-200	164	860	262	149	618	298	927	176	850	318	...	342	3859	234	134	258	2110	586	309	434	431
cons_length-250	164	859	265	297	618	298	1060	176	850	318	...	342	3860	234	134	258	2112	570	309	431	431
cons_length-300	164	859	264	297	618	495	1087	176	850	318	...	342	1350	234	134	258	2117	1671	309	433	431
cons_length-350	164	859	264	297	618	496	1087	176	850	331	...	342	4082	234	134	258	2118	1665	309	432	431
cons_length-400	164	859	264	297	618	496	3978	176	850	318	...	342	4075	234	129	258	2118	1665	309	432	431
cons_length-450	164	859	265	277	614	496	3967	176	854	319	...	342	4081	234	129	258	2116	1661	309	432	465
cons_length-500	164	859	265	277	614	496	3965	176	854	319	...	342	1360	234	194	258	2116	1662	309	432	465

9 rows × 40 columns

fig = px.line(cons_length)
fig.show()

The distance to the sample at 500 sequences

cons_div = data[['div-100', 'div-150', 'div-200', 'div-250', 'div-300', 'div-350', 'div-400', 'div-450', \
                    'div-500']]
cons_div = cons_div.T
cons_div

cluster	cluster_9	cluster_8	cluster_7	cluster_6	cluster_5	cluster_4	cluster_39	cluster_38	cluster_37	cluster_36	...	cluster_17	cluster_16	cluster_15	cluster_14	cluster_13	cluster_12	cluster_11	cluster_10	cluster_1	cluster_0
div-100	0.030120	0.000000	0.358634	0.220619	0.003247	0.397260	0.623491	0.369048	0.005869	0.081761	...	0.002924	0.509391	0.347826	0.262500	0.019231	0.010646	0.451852	0.247173	0.026713	0.328125
div-150	0.012195	0.001164	0.358634	0.227642	0.003247	0.399751	0.621832	0.011364	0.004695	0.080063	...	0.000000	0.508812	0.340517	0.310976	0.000000	0.008753	0.457330	0.000000	0.030093	0.325893
div-200	0.006098	0.001745	0.358634	0.305164	0.003247	0.413098	0.622240	0.011364	0.004695	0.076923	...	0.000000	0.508335	0.008547	0.310976	0.000000	0.006152	0.510676	0.000000	0.027714	0.328125
div-250	0.012195	0.001164	0.000000	0.062718	0.003247	0.410579	0.578905	0.011364	0.004695	0.070644	...	0.000000	0.509195	0.008547	0.310976	0.000000	0.006149	0.517921	0.000000	0.019699	0.325893
div-300	0.000000	0.001164	0.361059	0.062718	0.003247	0.015136	0.570467	0.011364	0.004695	0.070644	...	0.000000	0.008118	0.008547	0.310976	0.000000	0.004961	0.017102	0.000000	0.015029	0.325893
div-350	0.000000	0.000000	0.361059	0.062718	0.003247	0.004032	0.570071	0.005682	0.004695	0.255385	...	0.000000	0.522969	0.008547	0.310976	0.000000	0.003779	0.011121	0.000000	0.013889	0.325893
div-400	0.000000	0.000000	0.361059	0.059233	0.003247	0.004032	0.014478	0.005682	0.004695	0.004710	...	0.000000	0.521251	0.000000	0.287926	0.000000	0.003307	0.007514	0.000000	0.009259	0.330357
div-450	0.000000	0.000000	0.000000	0.003610	0.000000	0.004032	0.000504	0.005682	0.000000	0.003135	...	0.000000	0.522514	0.000000	0.294118	0.000000	0.001418	0.010533	0.000000	0.009259	0.000000
div-500	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	...	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000

9 rows × 40 columns

fig = px.line(cons_div)
fig.show()